(** Some needed definitions *)

module Inc = struct
  let some a = Some a
  let id a = a
  let id2 a b = a, b
  let id3 a b c = a, b, c
  let id4 a b c d = a, b, c, d
  let id5 a b c d e = a, b, c, d, e

  let maybe f = function Some a -> f a | None -> ()
  let maybe_opt f = function Some a -> Some (f a) | None -> None
end

open Inc
let scan_float txt = Scanf.sscanf txt "%f" id
let scan_int txt = Scanf.sscanf txt "%i" id
external string_to_format :
           string -> ('a, 'b, 'c, 'd, 'e, 'f) format6 = "%identity"

module Bin = struct
  let rec (--) i1 i2 = if (i1>i2) then [] else i1::((i1+1)--i2)
  let ( += ) acc c = acc:= !acc +. c
(*  let (|-) i1 i2 = Array.init (i2-i1+1) (fun i -> i+i1);;*)
  let (@*) a b = a, b
  let ($) a b x = a (b x)
(*  let ( // ) dir file = 
      ( match Sys.os_type with
        | "Unix" | "Cygwin" -> 
                dir ^ "/" ^ file
        | "Win32" -> 
                dir ^ "\\" ^ file
        | _ ->
                failwith "//: unknown os type"
      )
*)
  let ( // ) = Filename.concat
  let ( !// ) dir = 
      ( match Sys.os_type with
        | "Unix" | "Cygwin" -> 
                "/" ^ dir
        | "Win32" -> 
                "C:\\" ^ dir
        | _ ->
                failwith "!//: unknown os type"
      )
end

module Cycle = struct
  type ('a, 'b) result = [ `Result of 'b | `Continue of 'a ]
  let rec repeat ~from:(a0: 'a) (f: 'a -> ('a,'b) result) : 'b =
    match (f a0) with 
    | `Continue a1 -> repeat ~from:a1 f 
    | `Result x -> x
  ;;
  type 'a noresult = [ `End | `Continue of 'a ]
  let rec iter ~from:(a0: 'a) (f: 'a -> 'a noresult) =
    match (f a0) with 
    | `Continue a1 -> iter ~from:a1 f 
    | `End -> ()
end

let print10 name array =
    print_string name;
    Array.iteri (fun i a ->
        Printf.printf "\t%.2f " a;
        if i mod 10 = 9
        then print_newline ();
    ) array;
    print_newline ()

let iter64 ~range ~(f: Int64.t -> unit) =
    let first, last = range in
    let n = ref first in
    while !n < last do
        f !n;
        n := Int64.add !n 1L
    done
  ;;
let is_infinite x = classify_float x = FP_infinite

(* let is_nan x = classify_float x = FP_nan *)
let is_nan x = (x:float) <> x (* vastly faster; this is important not
                                 to slow down the usual case where the
                                 result is not NaN *)

(*let not_nan x = not (is_nan x)*)
let not_nan x = (x = x)

(*
let is_finite x =
  let x = classify_float x in
  x = FP_normal || x = FP_subnormal || x = FP_zero
*)

let is_finite x = (x <> infinity) && (x <> neg_infinity)
        && (not_nan x)

module Make_safe_float_op(
    Cont: sig
        type t 
        type place_info
        val fold: ('a -> float -> 'a) -> 'a -> t -> 'a
        val iterk: (place_info -> float -> unit) -> t -> unit
    end) =
struct
  (** ignoring all nans from initial data *)

  let sum arr =
    Cont.fold (
        fun acc c -> if not_nan c then acc+.c else acc
     ) 0. arr
  ;;
  (* same, finite *)
  let sumf arr =
    Cont.fold (
        fun acc c -> if not_nan c && is_finite c then acc+.c else acc
     ) 0. arr
(*
  (** find minimum of array, and maximum *)
  let minmax arr =
    let result = Cont.fold (fun result a ->
        if is_nan a 
        then result
        else (
            match result with
                Some (min, max) -> 
                    Some ((if a < min then a else min),
                          (if a > max then a else max))
              | None ->
                      Some (a, a)
        )
    ) None arr 
    in match result with Some (min, max) -> (min, max) | None -> (nan,nan)

  (** same as minmax, but only finite values are allowed *)
  let minfmaxf arr =
    let result = Cont.fold (fun result a ->
        if is_nan a || a = infinity || a = neg_infinity
        then result
        else (
            match result with
                Some (min, max) -> 
                    Some ((if a < min then a else min),
                          (if a > max then a else max))
              | None ->
                      Some (a, a)
        )
    ) None arr 
    in match result with Some (min, max) -> (min, max) | None -> (nan,nan)
*)
  let place_of_c pred cond arr =
    let c = ref None in
    Cont.iterk (fun p v ->
        if cond v then (
            match !c with
              | Some (_, acc) ->
                if pred v acc then c := Some (p,v)
              | None ->
                  c := Some (p,v)
        )
    ) arr;
    match !c with Some (p,_) -> p | _ -> raise Not_found
  ;;

  let place_of_max = place_of_c (>) not_nan
  let place_of_min = place_of_c (<) not_nan
  let place_of_maxf = place_of_c (>) (fun a -> not_nan a && is_finite a)
  let place_of_minf = place_of_c (<) (fun a -> not_nan a && is_finite a)

  let calc_minmax_c cond arr =
    let result = Cont.fold (fun result a ->
        if cond a
        then (
            match result with
                Some (min, max) -> 
                    Some ((if a < min then a else min),
                          (if a > max then a else max))
              | None ->
                      Some (a, a)
        )
        else result
    ) None arr in
    match result with Some (a,b) -> (a,b) | None -> raise Not_found
  ;; 

  let minmax = calc_minmax_c not_nan
  let minfmaxf = calc_minmax_c (fun a -> not_nan a && is_finite a)
  let min0max = calc_minmax_c (fun a -> not_nan a && a > 0.)
  let min0maxf = calc_minmax_c (fun a -> not_nan a && a > 0. && is_finite a)
  (*
(** find minimum among positive values of array, and maximum *)
let min0max arr =
    Array.fold_left (fun result a ->
        if is_nan a || a <= 0.
        then result
        else (
            match result with
                Some (min, max) -> 
                    Some ((if a < min then a else min),
                          (if a > max then a else max))
              | None ->
                      Some (a, a)
        )
    ) None arr 

(** same as min0maxf but only finite values allowed allowed *)
let min0maxf arr =
    Array.fold_left (fun result a ->
        if is_nan a || a <= 0. || a = infinity
        then result
        else (
            match result with
                Some (min, max) -> 
                    Some ((if a < min then a else min),
                          (if a > max then a else max))
              | None ->
                      Some (a, a)
        )
    ) None arr 
(*    let max=ref arr.(0) in
    for k=1 to Array.length arr - 1 do
        let a=arr.(k) in
        if a > !max then max := a;
    done;
    let min=ref !max in
    for k=1 to Array.length arr - 1 do
        let a=arr.(k) in
        if a > 0. && a < !min then min := a;
    done;
    !min, !max*)
*)
end

include Make_safe_float_op (struct
    type t = float array type place_info = int
    let fold = Array.fold_left let iterk = Array.iteri
end)

let gaus const mean sigma=
    function x -> let x_m=x-.mean in let sigma2=sigma*.sigma in
    const *. exp(-.x_m*.x_m/.sigma2)

